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How to calculate the exponential of matrices in an explicit manner is one of 
fundamental problems in almost all subjects in Science. 
§ : "^^P^""'^ " Mathematical Physics or Quantum Optics many problems are re- 

^^' duced to this calculation by making use of some approximations whether they are 

ij appropriate or not. However, it is in general not easy. 

3 . In this paper we give a very useful formula which is both elementary and getting 

on with computer. 
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To calculate the exponential of matrices 

°° Am 



m=0 

is one of fundamental problems in almost all subjects in Science ^. 

In fact, in Mathematical Physics or Quantum Optics many problems are reduced to 
this calculation by use of some approximations. However, this is a very hard problem. 
See, for example, hard and "painful" calculations in |2j and j3j. 

In usual textbooks of Linear Algebra (see for example |3]) it is recommended to diag- 
onalize A like 

A = UDaU-^ for some U G GL{n; C) (2) 

where Da is the diagonal matrix consisting of all eigenvalues of A. Unfortunately, this 
method is not realistic as known well. 

Let us explain in more detail. First of all we write the characteristic equation (poly- 
nomial) of A : 

= \XE-A\ = X^+ piA"-i + ■ ■ ■ + Pn-iX + pn (3) 

where pi = —tiA and pn = (— l)"det(y4). (jH)) can be decomposed into 

A" + piA"-^ + ■ ■ ■ + Pn-iX + Pn={X- ai) (A - aa) ■ • • (A - a„) (4) 

where aj G C. From this we have 



Pi 






P2 = {-ly^^aiOj 



i<j 



n 

Pn = (-irn«^- (5) 



i=i 



^See the dictionary [Q concerning several mathematical notations in the paper 



For the eigenvalue a^ (j = 1, 2, ■ ■ ■ , n) we define \aj) G C" to be the eigenvector which 
is not necessarily normalized ^ 

A\aj) = aj\aj), (6) 

then we obtain the diagonalization (j21) if we define U as 



t/ = (|ai), |a2),--- , l«n)) e GL{n;C) 

and Da = diag(ai, a2,- ■ ■ , an). 

Weak points of this method are that 

(a) to find the eigenvectors © explicitly, 

(b) to calculate U~^ the inverse of U. 

They become more and more difficult as n becomes large. 



(7) 



Now, it is better to change the strategy. The famous theorem of Cayley-Hamilton 



gives 



A" + piA'^-' + ■■■+ pn-iA + PnE = ^A^ = -piA""-' p„_iA - pnE (8 



, so we can formally write 



A 



^-^ ml 



(9) 



m=0 



where fj = fj{pi,P2,''' ,Pn) for simplicity. The purpose in the following is to 
determine {/o, /i, ■ ■ ■ , fn-i} explicitly ^. 
Here we use the notation like 



foE + f,A + f2A^ + --- + U,,A 



n-l 



Efo + Ah + AV2 + ■ ■ ■ + A^- Vn-i = {E, A,A',---, A--') 



( f \ 

Jo 

h 

/2 



■ (10) 



\ /n-l / 

^In general, to normalize a system of vectors is not easy. See the example in the latter half. 
•^This way of thinking is very natural because A is finite-dimensional 



Note that 



A"= {E,A,A^,--- ,A""i) 



-Pn-1 



-P2 



from (jHl). 

From (jH)) we define the matrix 

/ 





1 
1 



-Pn 
-Pn-1 
-Pn-2 



\ 



\ 



■■■ 
■■■ 1 



:ii) 



-P2 

'Pi ) 

which is called the companion matrix. It is well known that L also satisfies the same 
characteristic equation 



n-l 



0= \\E-L\ = A"+piA"~^ + ---+p„_iA + p„. 

We reconsider an important role that L plays. 

Now, we are in a position to state the main result. 

Fundamental Lemma For any m G N U {0} we have 



yl™= {E,A,/^,--- ,A"-i)L'"ei 



rp 

where ei = (1, 0, 0, ■ ■ • ,0) (see the notation in (|TU)) once more). 
The proof is by mathematical induction. 



(12) 



(13) 



This leads us to 



e^= {E,A,A^,--- ,A"-i)eV. 



(14) 



That is, the calculation of e"^ is reduced to that of e^. We must again calculate the 
exponential ! What is interesting ? In this case, we can make L diagonal completely 
because it is simple enough ^. 
Let us solve the equation(s) 



L\aj) = aj\aj) for j = 1 ~ n , 



(15) 



which is easily obtained to be 



\a\ 



Pn-2 + Pn-3a H h pia''^^ + a""^ 



V 



Pi + a 

1 



/ 



(16) 



where a = ai,a2,- ■ ■ ,an ior simplicity. It is of course {ai\aj) ^ %. Then we have 



f/i = (lai), |a2),--- , |a„)) G GL{n;C) 



(17) 



and 



L = UrDAUr-^ 



^ e^ = ULe^'^UL''. (18) 

At first sight, it seems difficult to calculate Ul^^. However, we have a simple decom- 
position of Ul like 

Ul = PlQl (19) 

where P^ and Ql are given respectively as 



*There is another diagonalization (communicated by T. Suzuki) 



Pl 



'n-l 


Pn-2 


Pn-3 




Pl 


i 


'n-2 


Pn-3 




Pl 


1 





'n-3 














Pl 













Pl 


1 













1 





. . . 












l\ 



(20) 



Qi 



Ctl 



at 



^2 



CKn 



a„_i a^ 



a. 



n-l 



at 



tti 



, n-2 ^n-2 

an 



a 



n-2 ^n-2 



n-l 



a: 



,71—1 „n~l 
^2 



n— 1 



(21) 



\ "l "2 • • • «n-l "n. / 

Therefore t/^^^ = Qj^^Pj^^. Must we calculate both Q^^ and P^^ again ? From (fT^ it is 
just U^^ei not [/^^ itself that we must calculate. It is easy to see 



where e„ = (0, 0, ■ ■ ■ , 0, 1) and Q2 e„ is given as 
' the cof actor of a"^^ ' 



Ql'e 



L ^n 



\Qi 



the cofactor of d^ ^ 

the cofactor of a^l} 
the cofactor of a""^ 



1) 



71+1 



/ 



1 



n-=i,,^i(".-i) 

1 

n"=i,j>£2K-«2) 



n-^i 



3^^,3i'-n.- 



j(qj-q„_i) 



Moreover, we can determine Ul completely. From (jSJ we define 



{pj)j. = Pj — {all the terms containing a^} for 1 < k < n. 



For example, 



(22) 



(23) 



(24) 



(Pi)i = -(a2H ha„), {pn-i)i = {-IT "2 



(jLyi 



(Pl)2 = -(«1 + "3 H ha„), (Pn- 



IJ2 



-1)" aia3---a„. 



Then 



/ 



Ul = (|ai), |«2),--- , |ttn)) 



(Pn-l)l 


(Pn-l)2 ■ 


■ iPn^l)n-l 


(Pn-l)„ 


bn-2)i 


(Pn-2)2 ■ 


■ {Pn^2)n-1 


(Pn-2)„ 


(Pl)l 


(Pl)2 ■ 


■ bl)n-l 


bl)n 


1 


1 


1 


1 



(25) 



Noting 



e^ei = f/ie^^f/i-^ei 



/ 



\n+l 



(Pn-l)l 


(Pn-l)2 ■ 


■ (Pn-l)n-l 


(Pn-l) 


(Pn-2)i 


(Pn-2)2 ■ 


■ iPn^2)n-l 


(Pn-2) 


bi)i 


(Pl)2 ■ 


■ (Pl)n-l 


(Pl)n 


1 


1 


1 


1 



we finally obtain 



n7=i.i^i{".-"i) 

e^^2 

n"=ij7^2("j-"2) 



n?=l.,-^„-l(aj-"n-l) 









(26) 



e-^ = [E, A,A\--- , A"-i) e^ei ~ foE + i\A + f^A^ + . . . + f^_,A 



n—l 



with 

in a complete manner ^. 



< / <ri- 1 



(27) 



That is, our calculation is based on only simple operations like the powers of matrices or 
sum of them, etc, which is easily performed by computer. 



Let us list some important examples (the case of n = 3 and 4). 



^The result is also obtained by the spectral decomposition of matrices (communicated by T. Suzuki). 
However, it seems to us that the method is much elementary. The result has not been written in standard 
textbooks in Linear Algebra as far as we know 



n = 3 For 



A 



I an ai2 ai3 \ 

fl21 0,22 fl23 
\ ^31 ^32 a33 y 



e M(3;C) 



(28) 



we have 



= \XE-A\ 

= X — (an + a22 + a33)A + (aiia22 + aiia33 + a22a33 — ai2a2i — 013031 — a23a32)A — 
(^11022033 + 012^23031 + 013O21O32 — O13O22O31 — 012021033 — 011023032). (29) 

The Cardano formula (see for example ■5j) gives three solutions {01,02,0:3}. Then 



Pl = — (Oi + O2 + O3), P2 = O1O2 + O1O3 + O2O3, P3 = — O1O2O3 



(30) 



and 



and 



and 



Ur 



//o^ 
/i 



e^ei 



L = 



^ -J93 \ 
1 -P2 

V 1 -Pl J 



= UlDaUl-' 



I P2+Piai + a\ p2+pia2 + al P2+Pi03 + o|\ 



Pl + Oi 

1 



Pi + a2 
1 



O2O3 O1O3 O1O2 

— (02 + 03) — (01 + 03) — (01 + 02 

1 1 1 



Pi + as 
1 

\ 



02O3 
(02 + 03) 
1 



01O3 
(oi + O3) 
1 



01O2 

(Oi + O2) 

1 



(31) 



(32) 





e°i 




(02- 


-ai)(a3- 

e°2 


-ai) 


(01- 


-a2)(a3- 

e°3 


-02) 



(ai— a3)(a2— as) 



(33) 



Therefore we have finally 



with 



/o 
/i 

/2 



foE + f,A + f,A' 



02036 



+ 



aia^e 



+ 



aia2e 



03 



(02 -ai)(a3 - oi) (ai - a2)(a3 - a2) (ai - a3)(a2 - 0:3) ' 
(a2 + a3)e"^ (ai + 03)6"^ («! + a2)e"=' 



(a2 - ai)(a3 - tti) (ai - a2)(a3 - "2) («! - a3)(a2 - a3) ' 



:,a2 



^"3 



+ 



+ 



(a2 -ai)(a3 -tti) (ai - a2)(a3 - "2) (ai - a3)(a2 - 03) ' 



(34) 



(35) 
(36) 
(37) 



n = A For 



A 



/ \ 

ail ^12 fll3 <2i4 

021 '3'22 Q.23 0-24 

^31 032 033 034 

fl41 042 O43 CL4.4. 



e M(4;C) 



we have 



(3J 



\XE-A\ 

X^ - (an + 022 + 033 + a44)A^ + (011022 + 011033 + 011044 + O22O33 

+0220-44 + O33O44 — O12O21 — O13O31 — O14O41 — O23O32 — O24O42 — 034043)A 

— (011022033 + O11O22O44 + O11O33O44 + O22O33O44 + O12O23O31 + O12O24O41 
+O13O21O32 + O13O34O41 + 014021042 + O14O31O43 + O23O34O42 + O24O32O43 

— 011023032 — O11O24O42 — O11O34O43 — O12O21O33 — O12O21O44 — O13O22O31 
~Oi303i044 — O14O22O41 — O14O33O41 — O22O34O43 — O23O32O44 — 024033042)A 

+ det{A) (39) 



where det(74) is omitted, see for example |3]. The Ferrari formula or Euler one (see j3]) 



gives four solutions {ai, 02, 03, a4} ^. Then 



and 



and 



Pl = — («! + 0^2 + «3 + a4), P2 = «ia2 + ai«3 + «ia4 + tt2«3 + «2«4 + a3tt4, 
]93 = — (aia2tt3 + CHtt2«4 + «ltt3tt4 + «2«3tt4), ^4 = aia2«3«4 (40) 



L 



f" 








-P4 


1 








-P3 





1 





-P2 


u 





1 


-Pl ) 



UlDaUl-' 



(41) 



' P3 + P2ai + pial + af ps + p2a2 + Pi^l + '^2 ^3+^203+^103 + 03 Pz + P2aA + Pi^l + a\ * 



P2 +PlCil + of 

Pl +ai 



P2 +PlCt2 + "2 
Pl + 02 



P2+Pl03 + a§ 

Pl + 03 



P2 +Pia4 + 04 
Pl +04 

1 



— a203a4 



-aia^a^ 



— aia204 



-Q;ia203 



02O3 + 0204 + 0304 aia3 + aia4 + 0304 Qia2 + oia;4 + 01204 aiQ2 + 0103 + 0203 

— (02 + 03 + 04) -(01 + 03 + 04) —(01 + 02 + 04) —(01 + 02 + 03) 

1111 



(42) 



and 



/ f \ 

Jo 

/i 

/2 
\f3 / 

Therefore we have finally 



c^ei 



-Ur 



e"i 



\ 



(o2-oi)(a3-ai)(a4-ai) 

e^ 

e^ 

{ai-a3){a2—a3){a4~a:i) 

e^ 

\ (Qi-04)(a2-a4)(a3-a4) / 



(43) 



foE + fiA + hA^ + hA' 



^They are of course too complicated 



(44) 



10 



with 



<^2«3tt4e"^ «iQ;3tt4e°^ 



/o = 7 ^/ ^. r + 



/l 



(a2 - «i)(a3 - ai)(a4 - "i) («! - tt2)(«3 - tt2)(«4 - 02) 

aia2tt4e"^ aia2a3e"* 

;i - as) (02 - tt3)(«4 - ^3) ("i - 04) ("2 - 04) ("3 - a 

[a2a3 + a2a4 + 0304)6"^ {ctids + aia4 + a3a4)e°'^ 



(02 — ai)(a3 — ai){a4 — ai) (ai — a2)(a3 — a2)(a4 — 02) 
(aia;2 + cna4 + 0204)6°'^ (010:2 + OiOs + a2a3)e"'* 



/2 



(ai - as) (02 - as)(«4 - as) ("i - 04) (^2 - "4) ("3 - "4) ' 
(a2 + "3 + a4)e"^ (ai + as + a4)e"2 



(46) 



(a2 - ai)(as - ai)(a4 - ai) (ai - a2)(a3 - a2)(a4 - a2) 
(ai + 0:2 + a4)e"^ (ai + 02 + a3)e"3 



(cKi — as)(a2 — as)(a4 — as) (ai — a4)(a2 — a4)(as — a4) 



/3 



(a2 — ai)(a3 — ai)(a4 — ai) (ai — a2)(a3 — a2)(a4 — a2) 



(ai - a3)(a2 - a3)(a4 - as) (ai - a4)(a2 - a4)(a3 - a4) 



(47) 



(4J 



A comment is in order. 

(1) At first sight, the eigenvalues in our formula ()44|) appear to be different all. However, 
it is not true. Let us explain this with a simple example. 
For A G M(2; C) 

0= \XE-A\ =A2 + j9iA + p2 = (A-a)(A-/3), 

which gives pi = —(a + /3), p2 = a/? and 




Then the formula is 

e^ = (£;,A)e^ 

In the following let us consider two cases. 
(I) ay^P: 

11 



1 




/3 _ 

In this case, L = Ul \ \ ^l with 

a 



. -a -(3 \ ^ 1/1/3 

^1 1 y /^-" \ -1 -a 

and it is easy to see 

(II) a = p: 

In this case, L is 

V 1 -Pi / V 1 2q 

we cannot diagonalize L, so we must use the method of Jordan canonical form. Since this 
method is well-known (see |4j) the details are omitted. If we set 



then it is easy to see 



Therefore 



-a 1 \ , / 1 

10/ Via 



a \ 1 2a 



a / \ e° 



we have 



, , 1 \ / (l-a)e^ 
eM J = I I • (50) 

Comparing (fSUj) with P^ it is clear that lim^q^^.^^ = (fSUj). That is, our method covers 
all special cases by taking some appropriate limit. 

12 



To obtain the formula including multi-plicities of eigenvalues is a good exercise. We leave 
it to readers. 

(2) In realistic problems we must calculate e*^ in place of e"^. The result is changed simply 
into 

e*^ = foE + f,A + f,A^ + ... + f^_,A^-^ (51) 

with 

^ ^ ( 1-j III p^^fc 

0</<n-l. (52) 



fi = { 
We leave the check to readers. 






(3) We comment on a generalization of the result. Let F be an entire function on C. 
Then the result is changed simply into 



F{A) = foE + fiA + f^A^ + ... + f^^.A 



n—l 



(53) 



with 



f=(_l)n+lSp {Pn-l-l)kF{ak) 



:t nj=i,jyfc(ai - "fc) 



< / <n- 1. 



(54) 



We leave the check to readers. 



We conclude this paper by stating our motivation. We are studying a quantum com- 
putation based on Cavity QED whose image is 





The general setting for a quantum computation based on Cavity QED 
the dotted line means a single photon inserted in the cavity and 
all curves mean external laser fields subjected to atoms 

13 



See 0; [ZI in detail. It is usually based on two level system of atoms. However, to take 
a multi level system of them into consideration may be better from the view point of 
decoherence which is a severe problem in quantum computation. To develop quantum 
circuits (see for example [H], [Hj and [Oj) we often encounter the problem to calculate the 
exponential of matrices in an explicit manner, which was very difficult (for at least Fujii). 
Since we have somewhat overcomed this difficulty it must be possible to reconsider 
quantum circuits in the multi level system. 
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